#ifndef MATH_3D_H
#define MATH_3D_H

#include "Prerequisites.h"

#include <math.h>
#include <vector>
#include <list>
#include <limits>


/** Wrapper class which indicates a given angle value is in Radians.
@remarks
	Radian values are interchangeable with Degree values, and conversions
	will be done automatically between them.
*/
class Degree;
class Vector2D;
class Vector3D;
class Vector4D;
class Quaternion;
class Ray;
class Plane;
class Sphere;
class AxisAlignedBox;
class Matrix4;


// interpolation functions
template<class T>
inline T interpolate(const float r, const T &v1, const T &v2)
{
	return static_cast<T>(v1*(1.0f - r) + v2*r);
}

template<class T>
inline T interpolateHermite(const float r, const T &v1, const T &v2, const T &in, const T &out)
{
	// basis functions
	float h1 = 2.0f*r*r*r - 3.0f*r*r + 1.0f;
	float h2 = -2.0f*r*r*r + 3.0f*r*r;
	float h3 = r*r*r - 2.0f*r*r + r;
	float h4 = r*r*r - r*r;

	// interpolation
	return static_cast<T>(v1*h1 + v2*h2 + in*h3 + out*h4);
}

enum Interpolations
{
	INTERPOLATION_NONE,
	INTERPOLATION_LINEAR,
	INTERPOLATION_HERMITE
};


class Radian
{
	float mRad;

public:
	explicit Radian ( float r=0 ) : mRad(r) {}
	Radian ( const Degree& d );
	Radian& operator = ( const float& f ) { mRad = f; return *this; }
	Radian& operator = ( const Radian& r ) { mRad = r.mRad; return *this; }
	Radian& operator = ( const Degree& d );

	float valueDegrees() const; // see bottom of this file
	float valueRadians() const { return mRad; }
	float valueAngleUnits() const;

	const Radian& operator + () const { return *this; }
	Radian operator + ( const Radian& r ) const { return Radian ( mRad + r.mRad ); }
	Radian operator + ( const Degree& d ) const;
	Radian& operator += ( const Radian& r ) { mRad += r.mRad; return *this; }
	Radian& operator += ( const Degree& d );
	Radian operator - () const { return Radian(-mRad); }
	Radian operator - ( const Radian& r ) const { return Radian ( mRad - r.mRad ); }
	Radian operator - ( const Degree& d ) const;
	Radian& operator -= ( const Radian& r ) { mRad -= r.mRad; return *this; }
	Radian& operator -= ( const Degree& d );
	Radian operator * ( float f ) const { return Radian ( mRad * f ); }
	Radian operator * ( const Radian& f ) const { return Radian ( mRad * f.mRad ); }
	Radian& operator *= ( float f ) { mRad *= f; return *this; }
	Radian operator / ( float f ) const { return Radian ( mRad / f ); }
	Radian& operator /= ( float f ) { mRad /= f; return *this; }

	bool operator <  ( const Radian& r ) const { return mRad <  r.mRad; }
	bool operator <= ( const Radian& r ) const { return mRad <= r.mRad; }
	bool operator == ( const Radian& r ) const { return mRad == r.mRad; }
	bool operator != ( const Radian& r ) const { return mRad != r.mRad; }
	bool operator >= ( const Radian& r ) const { return mRad >= r.mRad; }
	bool operator >  ( const Radian& r ) const { return mRad >  r.mRad; }
};

/** Wrapper class which indicates a given angle value is in Degrees.
@remarks
	Degree values are interchangeable with Radian values, and conversions
	will be done automatically between them.
*/
class Degree
{
	float mDeg; // if you get an error here - make sure to define/typedef 'float' first

public:
	explicit Degree ( float d=0 ) : mDeg(d) {}
	Degree ( const Radian& r ) : mDeg(r.valueDegrees()) {}
	Degree& operator = ( const float& f ) { mDeg = f; return *this; }
	Degree& operator = ( const Degree& d ) { mDeg = d.mDeg; return *this; }
	Degree& operator = ( const Radian& r ) { mDeg = r.valueDegrees(); return *this; }

	float valueDegrees() const { return mDeg; }
	float valueRadians() const; // see bottom of this file
	float valueAngleUnits() const;

	const Degree& operator + () const { return *this; }
	Degree operator + ( const Degree& d ) const { return Degree ( mDeg + d.mDeg ); }
	Degree operator + ( const Radian& r ) const { return Degree ( mDeg + r.valueDegrees() ); }
	Degree& operator += ( const Degree& d ) { mDeg += d.mDeg; return *this; }
	Degree& operator += ( const Radian& r ) { mDeg += r.valueDegrees(); return *this; }
	Degree operator - () const { return Degree(-mDeg); }
	Degree operator - ( const Degree& d ) const { return Degree ( mDeg - d.mDeg ); }
	Degree operator - ( const Radian& r ) const { return Degree ( mDeg - r.valueDegrees() ); }
	Degree& operator -= ( const Degree& d ) { mDeg -= d.mDeg; return *this; }
	Degree& operator -= ( const Radian& r ) { mDeg -= r.valueDegrees(); return *this; }
	Degree operator * ( float f ) const { return Degree ( mDeg * f ); }
	Degree operator * ( const Degree& f ) const { return Degree ( mDeg * f.mDeg ); }
	Degree& operator *= ( float f ) { mDeg *= f; return *this; }
	Degree operator / ( float f ) const { return Degree ( mDeg / f ); }
	Degree& operator /= ( float f ) { mDeg /= f; return *this; }

	bool operator <  ( const Degree& d ) const { return mDeg <  d.mDeg; }
	bool operator <= ( const Degree& d ) const { return mDeg <= d.mDeg; }
	bool operator == ( const Degree& d ) const { return mDeg == d.mDeg; }
	bool operator != ( const Degree& d ) const { return mDeg != d.mDeg; }
	bool operator >= ( const Degree& d ) const { return mDeg >= d.mDeg; }
	bool operator >  ( const Degree& d ) const { return mDeg >  d.mDeg; }
};

/** Wrapper class which identifies a value as the currently default angle
	type, as defined by Math::setAngleUnit.
@remarks
	Angle values will be automatically converted between radians and degrees,
	as appropriate.
*/
class Angle
{
	float mAngle;
public:
	explicit Angle ( float angle ) : mAngle(angle) {}
	operator Radian() const;
	operator Degree() const;
};

// these functions could not be defined within the class definition of class
// Radian because they required class Degree to be defined
inline Radian::Radian ( const Degree& d ) : mRad(d.valueRadians()) {
}
inline Radian& Radian::operator = ( const Degree& d ) {
	mRad = d.valueRadians(); return *this;
}
inline Radian Radian::operator + ( const Degree& d ) const {
	return Radian ( mRad + d.valueRadians() );
}
inline Radian& Radian::operator += ( const Degree& d ) {
	mRad += d.valueRadians();
	return *this;
}
inline Radian Radian::operator - ( const Degree& d ) const {
	return Radian ( mRad - d.valueRadians() );
}
inline Radian& Radian::operator -= ( const Degree& d ) {
	mRad -= d.valueRadians();
	return *this;
}

/** Class to provide access to common mathematical functions.
	@remarks
		Most of the maths functions are aliased versions of the C runtime
		library functions. They are aliased here to provide future
		optimisation opportunities, either from faster RTLs or custom
		math approximations.
	@note
		<br>This is based on MgcMath.h from
		<a href="http://www.geometrictools.com/">Wild Magic</a>.
*/
class Math
{
public:
   /** The angular units used by the API. This functionality is now deprecated in favor
	   of discreet angular unit types ( see Degree and Radian above ). The only place
	   this functionality is actually still used is when parsing files. Search for
	   usage of the Angle class for those instances
   */
   enum AngleUnit
   {
	   AU_DEGREE,
	   AU_RADIAN
   };

protected:
   // angle units used by the api
   static AngleUnit msAngleUnit;

	/// Size of the trig tables as determined by constructor.
	static int mTrigTableSize;

	/// Radian -> index factor value ( mTrigTableSize / 2 * PI )
	static float mTrigTableFactor;
	static float* mSinTable;
	static float* mTanTable;

	/** Private function to build trig tables.
	*/
	void buildTrigTables();

	static float SinTable (float fValue);
	static float TanTable (float fValue);
public:
	/** Default constructor.
		@param
			trigTableSize Optional parameter to set the size of the
			tables used to implement Sin, Cos, Tan
	*/
	Math(unsigned int trigTableSize = 4096);

	/** Default destructor.
	*/
	~Math();

	static inline int IAbs (int iValue) { return ( iValue >= 0 ? iValue : -iValue ); }
	static inline int ICeil (float fValue) { return int(ceil(fValue)); }
	static inline int IFloor (float fValue) { return int(floor(fValue)); }
	static int ISign (int iValue);

	static inline float Abs (float fValue) { return float(fabs(fValue)); }
	static inline Degree Abs (const Degree& dValue) { return Degree(fabs(dValue.valueDegrees())); }
	static inline Radian Abs (const Radian& rValue) { return Radian(fabs(rValue.valueRadians())); }
	static Radian ACos (float fValue);
	static Radian ASin (float fValue);
	static inline Radian ATan (float fValue) { return Radian(atan(fValue)); }
	static inline Radian ATan2 (float fY, float fX) { return Radian(atan2(fY,fX)); }
	static inline float Ceil (float fValue) { return float(ceil(fValue)); }

	/** Cosine function.
		@param
			fValue Angle in radians
		@param
			useTables If true, uses lookup tables rather than
			calculation - faster but less accurate.
	*/
	static inline float Cos (const Radian& fValue, bool useTables = false) {
		return (!useTables) ? float(cos(fValue.valueRadians())) : SinTable(fValue.valueRadians() + HALF_PI);
	}
	/** Cosine function.
		@param
			fValue Angle in radians
		@param
			useTables If true, uses lookup tables rather than
			calculation - faster but less accurate.
	*/
	static inline float Cos (float fValue, bool useTables = false) {
		return (!useTables) ? float(cos(fValue)) : SinTable(fValue + HALF_PI);
	}

	static inline float Exp (float fValue) { return float(exp(fValue)); }

	static inline float Floor (float fValue) { return float(floor(fValue)); }

	static inline float Log (float fValue) { return float(log(fValue)); }

	static inline float Pow (float fBase, float fExponent) { return float(pow(fBase,fExponent)); }

	static float Sign (float fValue);
	static inline Radian Sign ( const Radian& rValue )
	{
		return Radian(Sign(rValue.valueRadians()));
	}
	static inline Degree Sign ( const Degree& dValue )
	{
		return Degree(Sign(dValue.valueDegrees()));
	}

	/** Sine function.
		@param
			fValue Angle in radians
		@param
			useTables If true, uses lookup tables rather than
			calculation - faster but less accurate.
	*/
	static inline float Sin (const Radian& fValue, bool useTables = false) {
		return (!useTables) ? float(sin(fValue.valueRadians())) : SinTable(fValue.valueRadians());
	}
	/** Sine function.
		@param
			fValue Angle in radians
		@param
			useTables If true, uses lookup tables rather than
			calculation - faster but less accurate.
	*/
	static inline float Sin (float fValue, bool useTables = false) {
		return (!useTables) ? float(sin(fValue)) : SinTable(fValue);
	}

	static inline float Sqr (float fValue) { return fValue*fValue; }

	static inline float Sqrt (float fValue) { return float(sqrt(fValue)); }

	static inline Radian Sqrt (const Radian& fValue) { return Radian(sqrt(fValue.valueRadians())); }

	static inline Degree Sqrt (const Degree& fValue) { return Degree(sqrt(fValue.valueDegrees())); }

	/** Inverse square root i.e. 1 / Sqrt(x), good for vector
		normalisation.
	*/
	static float InvSqrt(float fValue);

	static float UnitRandom ();  // in [0,1]

	static float RangeRandom (float fLow, float fHigh);  // in [fLow,fHigh]

	static float SymmetricRandom ();  // in [-1,1]

	/** Tangent function.
		@param
			fValue Angle in radians
		@param
			useTables If true, uses lookup tables rather than
			calculation - faster but less accurate.
	*/
	static inline float Tan (const Radian& fValue, bool useTables = false) {
		return (!useTables) ? float(tan(fValue.valueRadians())) : TanTable(fValue.valueRadians());
	}
	/** Tangent function.
		@param
			fValue Angle in radians
		@param
			useTables If true, uses lookup tables rather than
			calculation - faster but less accurate.
	*/
	static inline float Tan (float fValue, bool useTables = false) {
		return (!useTables) ? float(tan(fValue)) : TanTable(fValue);
	}

	static inline float DegreesToRadians(float degrees) { return degrees * fDeg2Rad; }
	static inline float RadiansToDegrees(float radians) { return radians * fRad2Deg; }

   /** These functions used to set the assumed angle units (radians or degrees)
		expected when using the Angle type.
   @par
		You can set this directly after creating a new Root, and also before/after resource creation,
		depending on whether you want the change to affect resource files.
   */
   static void setAngleUnit(AngleUnit unit);
   /** Get the unit being used for angles. */
   static AngleUnit getAngleUnit(void);

   /** Convert from the current AngleUnit to radians. */
   static float AngleUnitsToRadians(float units);
   /** Convert from radians to the current AngleUnit . */
   static float RadiansToAngleUnits(float radians);
   /** Convert from the current AngleUnit to degrees. */
   static float AngleUnitsToDegrees(float units);
   /** Convert from degrees to the current AngleUnit. */
   static float DegreesToAngleUnits(float degrees);

   /** Checks whether a given point is inside a triangle, in a
		2-dimensional (Cartesian) space.
		@remarks
			The vertices of the triangle must be given in either
			trigonometrical (anticlockwise) or inverse trigonometrical
			(clockwise) order.
		@param
			p The point.
		@param
			a The triangle's first vertex.
		@param
			b The triangle's second vertex.
		@param
			c The triangle's third vertex.
		@returns
			If the point resides in the triangle, <b>true</b> is
			returned.
		@par
			If the point is outside the triangle, <b>false</b> is
			returned.
	*/
	static bool pointInTri2D(const Vector2D& p, const Vector2D& a,
		const Vector2D& b, const Vector2D& c);

   /** Checks whether a given 3D point is inside a triangle.
   @remarks
		The vertices of the triangle must be given in either
		trigonometrical (anticlockwise) or inverse trigonometrical
		(clockwise) order, and the point must be guaranteed to be in the
		same plane as the triangle
	@param
		p The point.
	@param
		a The triangle's first vertex.
	@param
		b The triangle's second vertex.
	@param
		c The triangle's third vertex.
	@param
		normal The triangle plane's normal (passed in rather than calculated
			on demand since the callermay already have it)
	@returns
		If the point resides in the triangle, <b>true</b> is
		returned.
	@par
		If the point is outside the triangle, <b>false</b> is
		returned.
	*/
	static bool pointInTri3D(const Vector3D& p, const Vector3D& a,
		const Vector3D& b, const Vector3D& c, const Vector3D& normal);
	/** Ray / plane intersection, returns boolean result and distance. */
	static std::pair<bool, float> intersects(const Ray& ray, const Plane& plane);

	/** Ray / sphere intersection, returns boolean result and distance. */
	static std::pair<bool, float> intersects(const Ray& ray, const Sphere& sphere,
		bool discardInside = true);

	/** Ray / box intersection, returns boolean result and distance. */
	static std::pair<bool, float> intersects(const Ray& ray, const AxisAlignedBox& box);

	/** Ray / box intersection, returns boolean result and two intersection distance.
	@param
		ray The ray.
	@param
		box The box.
	@param
		d1 A real pointer to retrieve the near intersection distance
			from the ray origin, maybe <b>null</b> which means don't care
			about the near intersection distance.
	@param
		d2 A real pointer to retrieve the far intersection distance
			from the ray origin, maybe <b>null</b> which means don't care
			about the far intersection distance.
	@returns
		If the ray is intersects the box, <b>true</b> is returned, and
		the near intersection distance is return by <i>d1</i>, the
		far intersection distance is return by <i>d2</i>. Guarantee
		<b>0</b> <= <i>d1</i> <= <i>d2</i>.
	@par
		If the ray isn't intersects the box, <b>false</b> is returned, and
		<i>d1</i> and <i>d2</i> is unmodified.
	*/
	static bool intersects(const Ray& ray, const AxisAlignedBox& box,
		float* d1, float* d2);

	/** Ray / triangle intersection, returns boolean result and distance.
	@param
		ray The ray.
	@param
		a The triangle's first vertex.
	@param
		b The triangle's second vertex.
	@param
		c The triangle's third vertex.
	@param
		normal The triangle plane's normal (passed in rather than calculated
			on demand since the callermay already have it), doesn't need
			normalised since we don't care.
	@param
		positiveSide Intersect with "positive side" of the triangle
	@param
		negativeSide Intersect with "negative side" of the triangle
	@returns
		If the ray is intersects the triangle, a pair of <b>true</b> and the
		distance between intersection point and ray origin returned.
	@par
		If the ray isn't intersects the triangle, a pair of <b>false</b> and
		<b>0</b> returned.
	*/
	static std::pair<bool, float> intersects(const Ray& ray, const Vector3D& a,
		const Vector3D& b, const Vector3D& c, const Vector3D& normal,
		bool positiveSide = true, bool negativeSide = true);

	/** Ray / triangle intersection, returns boolean result and distance.
	@param
		ray The ray.
	@param
		a The triangle's first vertex.
	@param
		b The triangle's second vertex.
	@param
		c The triangle's third vertex.
	@param
		positiveSide Intersect with "positive side" of the triangle
	@param
		negativeSide Intersect with "negative side" of the triangle
	@returns
		If the ray is intersects the triangle, a pair of <b>true</b> and the
		distance between intersection point and ray origin returned.
	@par
		If the ray isn't intersects the triangle, a pair of <b>false</b> and
		<b>0</b> returned.
	*/
	static std::pair<bool, float> intersects(const Ray& ray, const Vector3D& a,
		const Vector3D& b, const Vector3D& c,
		bool positiveSide = true, bool negativeSide = true);

	/** Sphere / box intersection test. */
	static bool intersects(const Sphere& sphere, const AxisAlignedBox& box);

	/** Plane / box intersection test. */
	static bool intersects(const Plane& plane, const AxisAlignedBox& box);

	/** Ray / convex plane list intersection test.
	@param ray The ray to test with
	@param plaeList List of planes which form a convex volume
	@param normalIsOutside Does the normal point outside the volume
	*/
	static std::pair<bool, float> intersects(
		const Ray& ray, const std::vector<Plane>& planeList,
		bool normalIsOutside);
	/** Ray / convex plane list intersection test.
	@param ray The ray to test with
	@param plaeList List of planes which form a convex volume
	@param normalIsOutside Does the normal point outside the volume
	*/
	static std::pair<bool, float> intersects(
		const Ray& ray, const std::list<Plane>& planeList,
		bool normalIsOutside);

	/** Sphere / plane intersection test.
	@remarks NB just do a plane.getDistance(sphere.getCenter()) for more detail!
	*/
	static bool intersects(const Sphere& sphere, const Plane& plane);

	/** Compare 2 reals, using tolerance for inaccuracies.
	*/
	static bool RealEqual(float a, float b,
		float tolerance = std::numeric_limits<float>::epsilon());

	/** Calculates the tangent space vector for a given set of positions / texture coords. */
	static Vector3D calculateTangentSpaceVector(
		const Vector3D& position1, const Vector3D& position2, const Vector3D& position3,
		float u1, float v1, float u2, float v2, float u3, float v3);

	/** Build a reflection matrix for the passed in plane. */
	static Matrix4 buildReflectionMatrix(const Plane& p);
	/** Calculate a face normal, including the w component which is the offset from the origin. */
	static Vector4D calculateFaceNormal(const Vector3D& v1, const Vector3D& v2, const Vector3D& v3);
	/** Calculate a face normal, no w-information. */
	static Vector3D calculateBasicFaceNormal(const Vector3D& v1, const Vector3D& v2, const Vector3D& v3);
	/** Calculate a face normal without normalize, including the w component which is the offset from the origin. */
	static Vector4D calculateFaceNormalWithoutNormalize(const Vector3D& v1, const Vector3D& v2, const Vector3D& v3);
	/** Calculate a face normal without normalize, no w-information. */
	static Vector3D calculateBasicFaceNormalWithoutNormalize(const Vector3D& v1, const Vector3D& v2, const Vector3D& v3);

	/** Generates a value based on the gaussian (normal) distribution function
		with the given offset and scale parameters.
	*/
	static float gaussianDistribution(float x, float offset = 0.0f, float scale = 1.0f);

	static const float POS_INFINITY;
	static const float NEG_INFINITY;
	static const float PI;
	static const float TWO_PI;
	static const float HALF_PI;
	static const float fDeg2Rad;
	static const float fRad2Deg;

	//added functions
	static float frand();
	static float randfloat(float lower, float upper);
	static int randint(int lower, int upper);

	template<class T>
	static T lifeRamp(float life, float mid, const T &a, const T &b, const T &c);

	static bool isZero(Vector3D& v);
	static float getDegreeBetween(Vector3D& v1, Vector3D& v2);
	static float getRadianBetween(Vector3D& v1, Vector3D& v2);

	inline static float uint16ToFloat(const uint16 value) { return (float)value/65536.f; }

};

// these functions must be defined down here, because they rely on the
// angle unit conversion functions in class Math:

inline float Radian::valueDegrees() const
{
	return Math::RadiansToDegrees ( mRad );
}

inline float Radian::valueAngleUnits() const
{
	return Math::RadiansToAngleUnits ( mRad );
}

inline float Degree::valueRadians() const
{
	return Math::DegreesToRadians ( mDeg );
}

inline float Degree::valueAngleUnits() const
{
	return Math::DegreesToAngleUnits ( mDeg );
}

inline Angle::operator Radian() const
{
	return Radian(Math::AngleUnitsToRadians(mAngle));
}

inline Angle::operator Degree() const
{
	return Degree(Math::AngleUnitsToDegrees(mAngle));
}

inline Radian operator * ( float a, const Radian& b )
{
	return Radian ( a * b.valueRadians() );
}

inline Radian operator / ( float a, const Radian& b )
{
	return Radian ( a / b.valueRadians() );
}

inline Degree operator * ( float a, const Degree& b )
{
	return Degree ( a * b.valueDegrees() );
}

inline Degree operator / ( float a, const Degree& b )
{
	return Degree ( a / b.valueDegrees() );
}

template<class T>
T Math::lifeRamp(float life, float mid, const T &a, const T &b, const T &c)
{
	if (life<=mid) return interpolate<T>(life / mid, a, b);
	else return interpolate<T>((life-mid) / (1.0f-mid), b, c);
}

#endif
